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Analysis  of  Laminar  and  Turbulent  Symmetric 
Blunt  Tralllng-Edge  Flows 


SUMMARY 


Contributions  are  made  to  interacting  viscous  flow  theory  for  symmetric 
separated  flow  past  blunt  airfoil  trailing  edges.  A  theoretical  and  compu¬ 
tational  technique  for  predicting  turbulent,  high  Reynolds  number,  separation 
and  near  wake  flow  is  provided.  A  detailed  assessment  of  a  viscous/inviscld 
iteration  and  an  Inverse  boundary  layer  solution  procedure  is  provided  through 
a  systematic  series  of  test  examples.  These  involve  subsonic  laminar  and  tur¬ 
bulent  flows  past  flat-plate  and  thickened  flat-plate  (blunt  trailing-edged) 
airfoils.  In  addition,  a  general  approach  is  outlined  for  incorporating  local 
strong-interaction  solutions  into  an  overall  viscous  solution  procedure  for 
isolated  airfoils  and  cascades. 
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INTRODUCTION 


An  important  problem  faced  by  the  designer  of  advanced  gas  turbine 
engines  is  the  prediction  of  high  Reynolds  number  viscous  flow  and,  in  par¬ 
ticular,  viscous  separation  phenomena,  in  compressor  and  turbine  blade  pas¬ 
sages.  Viscous  effects  control  exit  flow  angles,  aerodynamic  losses,  heat 
transfer  rates,  and  blade  stall,  and  hence,  must  be  accounted  for.  Clearly 
then,  an  analytical  procedure  for  predicting  viscous  flow  in  high-performance 
compressor  and  turbine  blade  rows  would  be  a  significant  contribution  to  the 
successful  design  process. 

The  overall  goal  of  the  present  research  program  is  to  employ  interacting 
viscous/invlscld  flow  concepts  to  develop  a  reliable  analytical/nuraer ical 
method  for  predicting  viscous  effects  in  subsonic  to  low  supersonic  cascade 
flows  at  high  Reynolds  numbers.  The  approach  to  be  followed  is  similar 
to  that  being  successfully  applied  to  high  Reynolds  number  external  aerodynamic 
flows  wherein  viscous/invlscld  interaction  concepts  are  used  to  construct  the 
full  flow  field  with  locally  relevant  component  flow  solvers  (Ref.  1).  Thus, 
for  the  case  of  high  Reynolds  number  flow  in  a  cascade,  the  full  flow  field 
is  broken  down  into  the  four  major  categories  to  represent  1)  the  inviscid 
flow;  2)  attached  boundary  layer  flow;  3)  wake  flow;  and  4)  locally  interact¬ 
ing  (possibly  separated)  flow  such  as  occurs  at  a  blunt  trailing  edge.  Con¬ 
struction  of  a  general  cascade  flow  solver  involves  first,  the  development  of 
component  flow  solvers,  and  second,  the  matching  of  these  component  solvers 
into  an  overall  computational  procedure  to  produce  a  reliable  and  efficient 
general  flow  solver.  Therefore,  the  utility  and  reliability  of  this  cascade 
solver  is  critically  dependent  on  the  level  of  development  of  its  component 
members.  The  first  three  components  which  deal  with  the  inviscid  flow,  the 
boundary  layer  flow,  and  the  wake  flow  are  well  understood  fluid  dynamically 
and  solution  methods  are  at  a  relatively  mature  state.  Work  on  the  fourth 
component,  the  local  interaction  solver,  has  been  progressing  at  UTRC  under 
contract  with  the  Naval  Air  Systems  Command.  In  the  present  phase  of  this 
effort,  solution  techniques  have  been  developed  and  demonstrated  for  symmetric, 
blunt,  trailing-edge  flows. 

The  present  effort  is  a  continuation  of  the  work  initiated  by  Werle  and 
Verdon  (Refs.  2  through  5),  in  which  a  theoretical  basis  for  addressing 
separated  trailing-edge  flows  was  established  using  asymptotic  triple— deck 
concepts.  A  computational  method  for  solving  the  triple-deck  equations  was 
first  demonstrated  for  supersonic  flow  (Refs.  2  and  3)  and  then  extended  to 
subsonic  flow  (Refs.  4  and  5)  using  Carter's  (Ref.  6)  Iterative  technique. 
Thereafter,  the  theoretical  concepts  were  generalized  to  the  finite  Reynolds 
number  regime  through  Introduction  of  a  special  separated  flow  coordinate  sys¬ 
tem  (Refs.  4  and  5).  Turbulence  effects  were  accommodated  through  introduction 
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of  a  turbulent  version  of  the  Levy-Lees  transformation  (Ref.  4).  The  goal  of 
the  current  effort  has  been  principally  to  refine,  demonstrate  and  verify  the 
theoretical  and  computational  methods  through  application  to  symmetric  hlgh- 
Reynolds-number  viscous  flows  past  blunt  airfoil  trailing  edges.  In  addition, 
attention  has  been  given  to  comprehending  the  "local"  viscous/inviscid  inter¬ 
action  solver  into  a  general  aerodynamic  formulation  in  which  the  local  strong 
viscous/inviscid  interaction  solution  is  coupled  to  a  global  inviscid  solu¬ 
tion  for  cascade  flow  fields. 

Results  are  presented  here  for  a  wide  range  of  trailing  edge  flows. 

First,  the  basic  finite  difference  algorithm  is  critically  assessed  through 
comparison  with  independent  triple-deck  solutions  for  laminar  flat-plate 
trailing  edge  flow.  Second,  the  interacting  boundary-layer  concept  and  solution 
algorithm  is  verified  through  a  systematic  Reynolds  number  study  showing  the 
necessary  approach  to  the  sharp  and  blunt  trailing-edge,  triple-deck,  solu¬ 
tions  as  Re  -*  Third,  the  extension  of  the  approach  to  turbulent  trailing 
edge  flows  is  verified  by  comparison  with  flat-plate  experimental  data;  and 
finally,  the  solution  procedure  is  demonstrated  for  symmetric,  turbulent, 
blunt  trailing-edge  flows. 
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FINITE  REYNOLDS  NUMBER  INTERACTION  MODEL 


General  Concepts 

For  flows  of  practical  Interest  in  either  external  or  internal 
aerodynamics  the  Reynolds  number  is  usually  sufficiently  high  so  that  the  flow 
past  an  airfoil  or  blade  can  be  divided  into  two  regions:  an  "inner"  dissipa¬ 
tive  region  consisting  of  the  boundary  layer  and  the  wake,  and  an  "outer" 
inviscid  region.  The  principal  interaction  between  the  viscous  and  inviscid 
regions  arises  from  the  displacement  thickness  effect  leading  to  a  thickened 
semi-infinite  equivalent  body  and  causing  changes  in  surface  pressures.  If 
the  interaction  is  "weak",  the  viscous  effect  on  the  pressure  is  small  (i.e., 

of  (Re-1/2) 

in  laminar  flow  or  of  Oiln  Re)  1  in  turbulent  flow  (Ref.  7)),  and 
the  complete  flow  problem  can  be  solved  in  a  hierarchial  manner  (Ref.  8).  The 
first  step  is  to  determine  the  inviscid  flow  past  the  airfoil.  The  resulting 
pressure  distribution  is  then  imposed  on  the  viscous  layer  calculation  to 
determine  the  displacement  due  to  viscous  dissipation  and  hence,  the  effective 
shape  of  the  body.  The  displacement  effect  on  the  inviscid  flow  is  then  cal¬ 
culated  to  provide  an  updated  surface  pressure  distribution. 

The  foregoing  method  for  calculating  the  interaction  between  the  viscous 
and  inviscid  parts  of  the  flow  is  based  on  a  direct  hierarchy  between  the 
viscous  and  inviscid  regions  which  applies  as  long  as  the  disturbances  to  the 
inviscid  flow  due  to  viscous  displacement  effects  are  small.  However,  the  flow 
over  an  airfoil  Involves  both  a  weak  overall  interaction  arising  from  standard 
displacement  effects  and  also  from  wake  curvature  effects,  and  local  "strong" 
interactions  arising,  for  example,  from  boundary  layer  separation,  shock 
boundary  layer  intL’raction,  or  the  rapid  flow  acceleration  immediately  aft  of 
the  trailing  edge.  In  such  situations,  viscous  displacements  cause  substantial 
changes  in  the  local  inviscid  pressure  field.  The  concept  of  an  inner  viscous 
region  and  an  outer  inviscid  region  still  applies  but  the  direct  hierarchial 
structure  of  the  flow  breaks  down.  in  a  strong  interaction  region  the  hier¬ 
archy  changes  from  direct  (i.e.,  pressure  determined  by  the  inviscid  flow)  to 
inverse  (pressure  determined  by  the  viscous  layer)  and  this  change  must  be 
accommodated  in  developing  a  complete  solution  (Ref.  8). 

The  ultimate  goal  of  the  present  research  program  is  to  develop  a  numeri¬ 
cal  solution  procedure  for  high  Reynolds  number  cascade  flow.  Our  intent  is 
to  employ  an  interacting  boundary  layer  approach  in  which  the  flow  in  the  in¬ 
viscid  region  is  potential  and  the  viscous  flow  is  governed  by  Prandtl's  boun¬ 
dary  layer  equations.  Numerical  solutions  of  the  full-potential  boundary  value 
problem  for  subsonic  (Ref.  9)  or  transonic  (Ref.  10)  cascades  and  direct  solu¬ 
tions  (i.e.,  pressure  prescribed)  for  the  boundary  layer  equations  (Ref.  11) 
have  been  developed  to  a  relatively  mature  state,  but  efficient  and  accurate 
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solutions  for  local  strong  viscous/lnviscid  interaction  regions  are  still 
under  development.  In  the  present  study  attention  has  been  focused  on  the 
trailing-edge  problem.  Here,  the  strong  Interaction  arises  from  the  abrupt 
change  in  the  slip  condition  that  the  boundary  layer  experiences  at  the  term¬ 
ination  of  the  airfoil  surface.  This  leads  to  a  singularity  in  the  classical 
boundary  layer  solution  and  subsequent  local  breakdown  of  a  weak- interaction 
solution  procedure.  An  Important  feature  of  the  present  formulation  is  the 
use  of  a  shear  layer  coordinate  (see  Ref.  4)  which  departs  from  the  airfoil 
surface  upstream  of  the  trailing  edge  and  "reattaches"  on  the  wake  centerline 
(Fig.  1).  After  making  the  usual  thin  viscous  layer  assumptions  in  directions 
along  and  normal  to  the  shear  layer  coordinate  and  applying  Prandtl's  transposi¬ 
tion  theorem,  the  classical  form  of  the  boundary  layer  equations  are  retained. 
The  present  report  deals  only  with  the  symmetric,  blunt,  trailing-edge  problem. 
Future  work  will  be  aimed  at  the  extension  of  present  concepts  and  solution 
procedures  to  compressor  and/or  turbine  blade  trailing-edge  flows,  which  are 
strongly  asymmetric,  and  thus  require  a  consideration  of  wake  curvature  effects. 
The  analysis  described  herein  has  been  developed  with  this  general  problem  in 
mind . 

*  *  *  * 

We  consider  high  Reynolds  number  (Re  =  u.^  L  /p^  )  flow,  with  neglig¬ 

ible  body  forces,  of  a  perfect  gas  with  constant  specific  heats  and  Prandtl 
number  past  a  two-dimensional,  symmetric  airfoil.  For  symmetric  flows  it  is 
sufficient  to  determine  a  solution  only  in  the  upper  half-plane  of  the  flow 
field.  In  the  following  discussion  flow  variables  and  spatial  coordinates  have 
been  made  dimensionless.  Lengths  have  been  scaled  with  respect  to  the  length 
of  the  airfoil  (L*) ,  density,  velocity,  and  viscosity  with  respect  to  their 
free-stream  values  (.r*,  ui  and  yS,  respectiv|ly) ,  pressure  with  respect  to 
twice  the  free-stream  dynamic  pressure  (p^  u*  ) ,  and  temperature  with  respect 
to  the  square  of  the  free-stream  speed  divided  by  the  specific  heat  at  con¬ 
stant  pressure  (u*^/c^).  Here  the  superscript  *  denotes  a  dimensional  quantity, 
and  the  subscript  «>  refers  to  conditions  at  infinity. 


Inviscid  Region 


The  inviscid  flow  is  assumed  to  be  isentropic  and  Irrotational .  The 
velocity  potential,  <1',  is  governed  by  the  equation 


=  Va>-V  (V«t)2/2 


(1) 


The  speed  of  sound  propagation.  A,  the  fluid  pressure,  p,  density,  p,  and 
temperature,  T,  are  determined  by  Bernoulli's  equation  and  the  Isentropic 
relations;  i.e.. 


(MpoA)^=  =  I  - 


(2) 
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where  M  is  the  Mach  number  and  Y  Is  the  specific  heat  ratio.  The  Inviscid  flow 
is  determined  as  a  solution  of  Eqs.  (1)  and  (2)  subject  to  conditions  of  flow 
tangency  at  the  displacement  surface,.^;  i.e., 

V4>-n  =0  r-ii 


where  n  denotes  the  unit  outward  normal  vector,  and  uniform  flow  in  the  far 
field;  i.e. , 

y4>-“ev  OS  R  — 00 

(4) 

where  e  is  a  unit  vector  directed  downstream  along  the  airfoil  chord  (i.e., 
the  x-axis)  and  ^  is  a  position  vector  emanating  from  a  point  on  the  airfoil. 
The  displacement  surface  is  defined  by  the  equation 

2)(s,n)  =  n -  S(s)+  t{s)  =  0  (5) 

where  s  denotes  the  shear  layer  coordinate,  n,  is  the  coordinate  normal  to 
the  shear  layer  coordinate,  6(s)  is  the  viscous  displacement  thickness, 
and  n  =  -t(s)  defines  the  airfoil  surface.  The  displacement  thickness  is 
measured  from  the  airfoil  surface  or  wake  centerline  in  a  direction  normal 
to  the  s-coordinate  curve  as  shown  in  figure  1. 

Solutions  for  interacting  flows  are  determined  by  an  iterative  procedure 
which  requires  successive  solutions  for  the  viscous  and  inviscid  regions.  If 
the  maximum  distance  between  the  displacement  and  shear-layer  coordinate 
surfaces,  j6(s)  -  t(s) (  ,  is  small,  updated  inviscid  predictions  can  be 

determined  efficiently  by  employing  a  perturbation  analysis  (c.f.  Ref.  12). 

For  this  purpose  the  velocity  potential  is  expanded  in  an  asymptotic  series; 
i.e. , 


0 +€<*>,  +  <?(  €2)  =  +  (?{e^) 

where  e  'v^do  -  tl^j^jj).  In  addition,  Taylor  series  expansions;  e.g., 

are  applied  to  refer  conditions  at  the  displacement  surface,^,  to  the  shear 
layer  coordinate  surface,^.  Equations  governing  the  zeroth  and  first-order 
potentials  are  obtained  by  substituting  the  series  expansions,  Eqs.  (6)  and 
(7),  into  the  governing  inviscid  equations,  Eqs.  (1)  through  (4),  equating 
terms  with  like  powers  of  c,  and  neglecting  terms  of  second  and  higher  order 
in  e. 
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The  zeroth  order  potential,  Is  a  solution  of  the  full-potential  equa¬ 
tion  (Eqs.  (1)  and  (2)  with  $  and  A  replaced  by  4>y  and  A^,  respectively) 
subject  to  the  condition  of  flow  tangency  at  the  shear  layer  coordinate 
surface;  i.e.. 


n 


=  0 


(8) 


and  uniform  conditions  in  the  far  field 

7*0 -ex,  R-od 


(9) 


The  first-order  potential  is  a  solution  of  the  linear,  variable-coefficient, 
differential  equation 


=[V«1>0'V+ 


(10) 


where  the  coefficients  depend  on  the  zeroth-order  solution.  First-order 
conditions  on  the  shear  layer  coordinate  surface  and  in  the  far  field  are 
given  by 


-  n  =  (v<l>o  ?)^(8-t)  -  (8-t)-^ 


^0  ,  on# 


dx\ 


(11) 


and 


V</>  =  0,  R  —  00 


(12) 


where  t  is  a  unit  vector  tangent  toS  and  directed  downstream.  The  inviscid 
pressure  at  the  displacement  surface  is  given  by 

p|2j  =  [Po+«P, -*•^5-1)-^]^+ 

(13) 
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where 


Po  (M<^Ao) 


2y/(y-l) 


y/{y-  I) 


(14) 


and 

(ep,  +  (S-t)^°)^  -  V4>  +  (S-t)  {V<I>o)^/2]| 

^  (15) 

Thus  the  inviscid  pressure  at  the  displacement  surface,^,  is  determined  in 
terms  of  information  supplied  at  the  shear- layer  coordinate  surface, 

Realistic  turbomachinir>  flow  configurations  frequently  consist  of 
thick,  highly-cambered  blades  and  large  flow  deflections  due  to  high  blade 
loading.  Thus  the  general  inviscid  model  derived  above  will  be  required  for 
compressor  or  turbine  interaction  analyses.  Numerical  solutions  of  the  non¬ 
linear,  zeroth-order  (Ref,  9)  and  the  linear,  first-order  (Ref.  13)  potential 
equations  have  already  been  developed  for  cascade  applications.  These 
solutions  and  their  associated  computer  codes  can  be  readily  modified  and 
extended  for  vlscous/inviscid  interaction  studies.  However,  in  the  present 
effort  our  primary  concern  has  been  the  application  of  an  inverse  viscous 
solution  procedure  and  a  viscous/inviscid  iteration  procedure  for  calculating 
the  flow  in  local  strong-interaction  regions.  For  this  purpose  the  inviscid 
model  can  be  simplified  further  by  restricting  consideration  to  symmetric 
flow  past  a  thin  airfoil.  In  this  case  the  zeroth-and  first-order  boundary 
valui’  problems  can  both  be  approximated  by  the  cl;ssical  small-disturbance 
equations,  and  the  inviscid  pressure  distribution  at  the  displacement  sur¬ 
face  can  be  determined  by  a  Cauchy  integral  evaluation;  i.e.. 


p^^ls) 


(Po+eP|)^ 

“  77-  Lj_  S-;  s-c 

h'{C)+8'{C)  , 

/_,  -  s-r~ 


(16) 
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where  h(s)  is  the  distance  from  the  airfoil  chord  to  the  airfoil  surface  in 
the  direction  normal  to  the  shear-layer  coordinate,  and  the  prime  denotes 
differentiation.  Note  that  to  within  the  first-order  approximation 


P^(S)  =  P2,(X)  =  p 


/•«>  h‘(C)  +  S‘(C) 

-n-  •'-oc  x-C 


(17) 


I 


J 


where  the  airfoil  thickness,  h,  and  viscous  displacement  thickness,  6, 
distributions  are  now  measured  normal  to  the  airfoil  chord  line. 


Viscous  Layer 

The  flow  in  the  inner  or  viscous  region  is  assumed  to  be  governed  by 
Prandtl's  boundary  layer  equations  referred  to  the  shear— layer  coordinate 
direction.  A  detailed  derivation  of  these  equations  can  be  found  in  Ref.  2. 
After  applying  Prandtl's  transposition;  i.e.  ,  introducing  the  scaled  normal 
coordinate,  iV,  and  velocity  component,  v,  where 

n  =  Re'''^n+t)  . 

V  =  Re'^(v  +  u“)  (18) 


and  u  and  v  are  the  fluid  velocity  components  in  the  s  and  n-directions, 
respectively,  the  governing  equations  become 


dipu)  dipV) 
ds  d?\ 


--0 


(19) 


s -momentum 
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2 

Here  H=T+u/2  is  the  total  enthalpy,  the  viscosity,  y,  is  assumed  to  be  a 
function  of  temperature  alone,  a  is  the  Prandtl  number,  the  subscript  e  denotes 
properties  at  the  edge  of  the  boundary  layer,  and  the  effective  turbulent 
viscosity  terms  P'p  and  y^  are  defined  by 


H-t  rye 

+  o’/ie/o-j 


where  e  is  the  turbulent  eddy  viscosity,  is  the  turbulent  Prandtl  number, 
and  is  the  longitudinal  intermittency  factor. 

The  foregoing  equations  are  solved  subject  to  the  following  boundary 
conditions : 


Edge  Conditions  (n  -*■  “) 

u— Ue  ,H-He 

Airfoil  Surface  Conditions  (n  =  0,  s  £ 

u  =  V  =0 

and  either 

<3H  -1/2 

H  :  Hw(s)  or  =  Re  '  qy^,ls) 


(23) 


(24) 


(25) 


where  and  9^^(s)  are  prescribed  functions 

Wake  Centerline  Conditions  (n  =  0,  s  >  s^j.) 


<?U  _ 

an  ^  an 


(26) 


Finally,  the  displacement  thickness  of  the  viscous  layer  is  defined  by  the 
relation 

8(s)  /"(i-  -^)dn 

(27) 
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Turbulence  Model 


The  turbulence  model  employed  in  the  present  study  is  similar  to  the 
Cebeci-Smith  (Ref.  14)  eddy  viscosity  formulation  as  modified  for  flat-plate 
near-wake  flows  by  Cebeci  et  al.  (Ref.  15).  Thus  for  wall  boundary  layers 
the  eddy  viscosity  formulation  considers  inner  and  outer  turbulence  layers 
with 


V') 

2  1/2  ~  -ll  1  ~  ~ 

e  =  6:  =  0.16  D  Re'  pnu  '  —  ,  n<n, 

‘  '  dn  ' 

(28a) 

■V. 

< 

e  =  Bo  =  0.0168 /3ue^’'Re'''^  ^"*(1- )dn,  n>n| 

(28b) 

jii 

f\j  <\j 

and  the  two  forms  being  equal  at  n  =  n^.  Here 

D=  [i-expC-n/o)] 

i 

(29a) 

j 

where  the  length  scale,  a,  is  given  by 

> 

a  =  26^Re''''^[/>{l-ll.8p^)uT]‘' 

(29b) 

i 

and  the  pressure  gradient,  p'*’,  and  shear  velocity,  u^,  are  defined  by 

1 

i.  1 

(29c) 

f . . 

1  ^ 

t- 

-i/4r  du  /  i‘/2 

UT=Re  [/i.  /^Jtiro 

(29d) 

In  the  wake,  the  inner  region  is  further  subdivided  into  two  parts. 


r 

L 


For  0  <  n  <  n^ 


-1/2  -I  ~ 

ew  =  0.4  Re  /JUrp  fj- 


where  n^,  is  the  value  of  ft  such  that 


•^CQIC  =  l•lUTp[2.5ln(Re'^%TEM' n)  +  5.2] 
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where  is  the  computed  velocity  profile  in  the  wake  and  is  the  shear 

velocity  at  a  cusped  trailing  edge  or  at  some  specified  location  upstream 
of  the  trailing  edge  for  a  blunt  trailing  edge.  When  n2  <  n  <  n^,  equation 
(28a)  is  applied  with  D  =  1.  To  prevent  e^  from  exceeding  e^  at  large  dis¬ 
tances  downstream  of  the  trailing  edge,  once  Cq  =  e^  the  eddy  viscosity  for 
larger  values  of  s  is  taken  to  be  constant  across  the  wake. 


The  longitudinal  intermittency  factor,  y  ,  is  taken  to  be  the  same  as  that 
employed  by  Harris  (Ref.  16);  i.e., 

Yj  =  l-exp[-0.412{s-Sti)^/Xy]  (31a) 

where  s^.^  is  the  location  of  transition  initiation,  and  the  scale  length  of 
the  transition  process  is  given  by 


^y-r^o.rs  ^y^:0.25 


(31b) 


Vlscid/Inviscld  Iteration 


The  complete  flow  field  is  determined  as  a  solution  of  the  inviscid 
equations  subject  to  the  condition  of  flow  tangency  at  the  displacement  sur¬ 
face,  2>,  and  a  solution  of  the  boundary  layer  equations  which  approaches 
Inviscid  conditions  on  at  the  outer  edge  of  the  viscous  layer  (i.e.,  u  u^ 
and  H  -*•  H  as  n  -*•«>) .  Since  the  displacement  surface  must  be  determined  as 
part  of  the  overall  solution,  an  iterative  solution  procedure  is  required. 

It  has  been  found  that  a  semi-inverse  calculation  using  underrelaxation 
provides  a  relatively  efficient  method  for  determining  the  flow  in  strong 
viscous/inviscid  interaction  regions  (Ref.  6).  In  this  approach  the  inviscid 
and  viscous  equations  are  solved  for  a  prescribed  nth  displacement  thickness 
distribution,  6^  (s),  to  determine  the  (n  +  l)th  inviscid  pressure  distribution, 
p7'*’^(s),  at  the  displacement  surface  and  the  (n  +  l)th  pressure  distribution, 
pjjT^c(s),  in  the  viscous  layer.  The  (n  +  l)th  estimate  for  the  displacement 
thickness  is  then  obtained  by  underrelaxation  i.e.. 


8 


n+i 


(s)  =  8"(s) 


(32) 


where  w<i  is  the  relaxation  parameter.  The  process  is  continued  until 


8""''(s)-8"(s)||  <  €g 


(33) 


12 


R81-914986-5 


where  ||  {|  Is  a  suitably  chosen  norm;  e.g.. 


Here  s  =  and  s  =  s^  define  the  upstream  and  downstream  boundaries  of  the 
strong  interaction  solution  domain. 

The  foregoing  procedure  consists  of  direct  solutions  of  the  Inviscid 
equations  and  Inverse  solutions  of  the  viscous  layer  equations  (l.e.,  6  is 
prescribed  to  determine  both  the  inviscid  and  viscous  pressure  fields).  An 
important  advantage  in  applying  an  inverse  procedure  to  determine  the  viscous 
flow  field  is  that  the  boundary  layer  equations  can  be  solved  through  a  separa¬ 
tion  point  and  into  a  region  of  reversed  flow,  whereas  a  direct  (l.e.,  dp/ds  = 

,  Ug  dug/ds  is  prescribed  to  determine  6(s))  boundary  layer  solution  will 
break  down  at  a  separation  point. 


13 
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VISCOUS  REGION  SOLUTION  PROCEDURE 


In  the  present  study  emphasis  has  been  placed  on  the  development  of 
a  viscous/inviscid  iteration  procedure  and  an  inverse  viscous  solution 
procedure  for  f inite-Reynolds- number ,  compressible,  trailing— edge  flow. 

For  this  purpose  consideration  was  restricted  to  adiabatic  flow  at  unit 
Prandtl  numbers  (o=o-p=l)  past  a  symmetric,  thin,  blunt-edged  airfoil.  The 
iteration  procedure  between  the  viscous  and  inviscid  solutions  was  de¬ 
scribed  in  the  previous  section  and  the  inverse  viscous  layer  procedure  is 
presented  here.  The  general  viscous  layer  equations  are  first  expressed 
in  "turbulent"  Levy-Lees  type  variables  and  then  after  appropriate  simpli¬ 
fications  are  introduced,  the  viscous  solution  procedure  is  described  for 
the  special  case  of  adiabatic  flow.  The  computer  code  being  developed  in 
conjunction  with  the  present  analysis  has  been  structured  such  that  ex¬ 
tensions  to  general  aerodynamic  profiles  and  heat  conducting  flows  can  be 
readily  Incorporated  in  future  studies. 


"Turbulent".  Levy-Lees  Transformation 

The  viscous  layer  equations  (Eqs.  (19),  (20)  and  (21))  can  be  recast  in 
terms  of  similarity  type  variables  by  employing  a  turbulent  generalization 
of  the  Levy-Lees  transformation  as  introduced  by  Werle  and  Verdon  (Ref.  4). 
This  transformation  serves  to  minimize  the  growth  of  the  viscous  layer  in 
the  computational  domain  and  eliminate  density  from  the  formal  governing  equa¬ 
tions  by  introducing  new  longitudinal  and  normal  coordinates  which  are  func¬ 
tions  of  density  and  viscosity.  The  new  independent  variables  are  defined  by 


and 


C,  +  /  /^e^e^eg  ds 


(35a) 


Ue 

1=-^  f  pan 


(35b) 


Hence  the  derivatives  become 
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±  . 

(3s  ■  Pe  Me  Ue  9 


J-  +  12L  -L 

dC  ds  d7j 


(36a) 


and 


(36b) 


Here,  the  definition  of  C  is  a  simple  generalization  of  the  usual 
laminar  Levy-Lees  streamwise  variable  because  of  the  introduction  of 
the  function  g(s).  This  function  can  be  specified  (see  below)  to  help 
control  the  thickness  of  the  turbulent  layer  in  the  -plane. 


With  new  dependent  variables,  F(^,n),  G(C»n),  and  V(5,n)»  defined  by 
F  =  u/ug 
G  =  H/Hg 


(37) 


ds  J/  ds 

the  governing  equations  reduce  to  the  familiar  form 
continuity 

C  -  moment UF 

d  L  dF' 


(38) 


2CF  +  V  +  i9(F^-G)-  =  0 

di  drj  d-q  \  d-q! 


energy 


(39) 


(40) 


di  '  drj  d7j\. . drj^  d-q 

where  the  following  notation  has  been  introduced 


u|/  [Te(i+u|/Tg)] 

11  ^  (,+  iiL) 

ue  dC  2Te  ' 

(41) 
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Boundary  conditions  for  the  foregoing  equations  are  as  follows  (c.f.  Eqs. 
(23)  through  (26)): 

Edge  Conditions  ( n->-“) 


FiG  — I 


Airfoil  Surface  Conditions  (C<C„„.n=0) 
—  “■  ~  ■  TE’ 

F  =  V=  0 


(42) 


(43) 


and  either 

G  =  Hw(^)/He  or  ||  =  (2f/Re)‘"%^(f) 

where  H^,(C)  and  qy(C)  are  prescribed  functions,  and  p 
G  and  the  fluid  properties  at  the  edge  of  the  boundary  layer  as  shown 
below 


/(pueHe) 

is  a  function  of  F, 


Wake  Centerline  Conditions 


dF  .  dG 


drj  drj 


=  V 


=  0 


(45) 


Fluid  properties  at  the  edge  of  the  boundary  layer  (i.e.,  inviscid 
properties  at  the  displacement  surface,  can  be  expressed  in  terms  of 
the  free-streara  and  edge  Mach  numbers.  It  follows  from  Bernoulli's 
equation,  Eq.  (2),  that 


(MooUe/’Ma)  =  (MooAe)* 


(46) 


and  therefore 


2|_  ^ 

^  ■  Me  dC  (^7b) 


In  addition,  it  follows  from  the  definition  of  displacement  thickness,  Eq. 
(27),  the  equation  of  state  (p=(y-1)pT/y) ,  and  the  boundary  layer  approx¬ 
imation  for  normal  pressure  gradient  (i.e.  3p/9n'V'0),  that 
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m=  PeUe8(Re/20'^*  =  ^"[6-f  +  !^  m|(G-f2)]  d-,^ 


where  m  Is  a  mass  flow  parameter  and 


y  -I  2 
Pelp  =(3  +  ^M, 


[g-f*] 


Thus  the  explicit  dependence  on  density  can  be  removed  from  the  surface 
heat  transfer  condition,  Eq.  (44).  Finally  a  shear  stress  parameter,  Cj 
is  defined  as  follows 

When  q  (or  n)  is  normal  to  the  airfoil  surface,  S,  Cf  represents  the 
surface  skin  friction  coefficient. 

Turbulence  Model 

The  eddy  viscosity  model,  outlined  in  the  previous  section  (c.f.  Eqs. 
(28)  through  (30),  is  rewritten  in  terms  of  the  turbulent  Levy-Lees  vari- 
able§  as  follows.  For  wall  boundary  layers 

e  =  Bj  =  0.16  D^(2CRe)'^^/5^  I  ,  ?)<  (5ia) 

e  =eo  =  0.0168 d?),  17 >7), 

where  is  the  value  of  q  such  the  ej^(C,q)  =  ejj(5,q)  and 


e  =  jf  d7]/p. 


D  =  i-exp  - 


a  =  26;i./[/3Re'^^(i-ii.8p^)uJ  ?■*■=  {MUe)^g  /(Reu^l 


Ut=  Ue{2€Re)''^‘*[^-~]^''J^  = 


For  wakes 


/  /  .►  \\/Z  rOoFlCiO)  T 

e  =e^  =  0.4(^ue)  (2fRe)  Ut^  exp[  -  5.2/2.5J  . 


(53a) 


J 


R81-914986-5 


i72<‘7<77, 


e  =eo. 


(53b) 

(53c) 


where  n2  ai'd  are  the  values  of  n  such  that  e^=e^  and  e^=e^,  respec¬ 
tively,  and  the  eddy  viscosity  is  taken  as  e^  across  the  wake  if  e^>eQ. 
Note  that  Eq.  (53a)  represents  an  approximation  to  the  model  of  Eqs. 

(30a)  and  (30b)  which  is  obtained  by  assuming  that  flow  properites  are 
constant  and  equal  to  their  centerline  values  over  the  inner  wake  region. 

The  transformation  function  g  is  chosen  so  that  the  coefficient  I 
defined  in  Eq.  (41)  is  approximately  equal  to  one  over  most  of  the  viscous 
layer.  This  is  accomplished  for  example  by  setting 

9"'+7t®o  (54) 

Thus  C.=pi;/ (ppUp)  for  both  laminar  flow  (y.j,  =  0)  and  in  the  outer  region  of 
a  fully-developed  turbulent  layer. 


Inverse  Solution  Procedure 

For  the  purpose  of  establishing  an  inverse  viscous  solution  procedure, 
the  following  simplifications  can  be  introduced  while  still  retaining  the 
essential  features  of  the  viscous  boundary  value  problem.  The  flow  is 
assumed  to  be  adiabatic  and  the  Prandtl  numbers  a  and  Oj  are  assumed  to  be 
unity.  Thus  G((l,n)  =  1,  and  it  is  only  necessary  to  solve  the  continuity 
and  momentum  equations  to  determine  the  functions  F  and  V.  In  addition,  the 
viscosity,  u*,  is  assumed  to  vary  linearly  with  temperature,  T  ,  in  the 
viscous  layer  and  therefor?. 


(55) 

and 

p/xU-t-y-e) 

.  «  =  - - -  (56) 

The  complete  flow  field  is  determined  by  an  iterative  procedure  (see 
Ref.  6)  which  is  continued  until  the  inviscid  solution  at  the  displacement 
surface,  ^  ,  matches  the  viscous  solution  at  the  edge  of  the  viscous  layer. 
At  each  step  of  this  process  the  viscous  layer  equations  are  solved  numer¬ 
ically  for  prescribed  mass  flow  parameter  (or  displacement  thickness;  c.f. 
Eq.  (48))  and  edge  Mach  number  distributions.  The  continuity.  Eq.  (38), 
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and  momentum  (Eq.  (39)  with  G  =  1)  equations  are  replaced  by  a  set  of  linear 
algebraic  equations  using  a  finite  difference  approximation  in  which  non¬ 
linear  terms  in  the  momentum  equation  are  quasi-linearized,  and  C  and  n-deriv- 
atives  are  replaced  by  one-sided  and  central  difference  expressions,  respec¬ 
tively  (see  Ref.  4).  Solutions  to  the  linear  difference  equations  are  deter¬ 
mined  by  a  superposition  technique;  i.e.,  new  dependent  variables 


F  =  +  ^TL 

V  =  /SVi  +  Vji 


(57) 


where 


Fj  =  0  and  Fjj  M  ,  17  —  00 

Fj  =  Vj  =  0  and  Fjj  =  Vjj  =  0  ,  £  <  C ,  T)  =  0 


-T —  and  —z —  =  0 , 
dr]  dr]  ' 


(58) 


are  introduced  so  that  explicit  dependence  on  the  pressure  gradient  par¬ 


ameter,  K,  is  eliminated  from  the  systems  of  equations,  I  and  II,  govern¬ 


ing  the  dependent  variables  Fj  and  Vj,  and  Fjj  and  V^j,  respectively. 
The  linear  systems  I  and  II,  are  each  solved  using  the  Davis  inversion 
scheme  (Ref.  17). 


To  recover  the  full  viscous  layer  solution,  6  must  be  determined  for 
the  given  mass  flow  parameter,  m'^(C),  and  edge  Mach  number,  M”(5),  dis¬ 
tributions  where  the  superscript  denotes  the  nth  iterative  estimate.  It 
follows  from  Eqs.  (48  with  G  =  1)  and  (57)  and  some  algebra  that  3  satis¬ 
fies  the  quadratic  equation 
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After  eliminating  the  extraneous  root,  the  solution  of  Eq.  (60)  is 


/3  =  (Ib/2Ia)[[l  +  4iQ(lc-m'^/l^]'/2-|] 


(61) 


Once  3  is  determined  the  (n  +  l)th  estimate  (viscous)  of  the  fluid  pro¬ 
perties  at  the  edge  of  the  viscous  layer  are  obtained  from  the  relation 
(c.f.  Eq.  (47b)) 


M* 


n-n 


(e)  =  MeU,)exp[  dc] 


2$ 


(62) 


where  the  subscript  1  denotes  the  initial  station  of  the  computational 
region, and  Bernoulli's  equation,  Eq.  (46).  The  pressure  determined  in 
this  manner:  i.e.,  p^l"^  ,  is  used  in  conjunction  with  the  inviscid  pres- 
sure  determined  by  Cauchy  integral  evaluation;  i.e.,  Eq.  (17)  with 
6(x)  =  6'^(x),  to  obtain  an  updated  displacement  thickness  distribution 
(c.f.  Eq.  (32)). 
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NUMERICAL  ALGORITHM 


The  basic  numerical  algorithm  for  solving  the  equations  described  in  the 
previous  sections  of  this  report  has  been  presented  in  detail  by  Werle  and 
Verdon  (Ref.  4).  Here  we  seek  mainly  to  identify  those  modifications  and/or 
additions  to  that  algorithm  required  for  treating  turbulent,  separated,  trailing- 
edge  flows.  First,  the  Cauchy  solver  for  the  inviscid  pressure  distribution  is 
modified  to  accommodate  a  variable  grid  spacing  in  the  streamwise  direction. 
Second,  modifications  to  the  difference  equation  coefficients  for  the  viscous 
equations  are  identified  for  accommodating  variable  grid  spacing  in  both  the 
streamwise  and  normal  directions.  Finally,  the  details  of  the  overall  solution 
algorithm  including  the  mesh  size  and  extent  considerations,  initialization 
procedures,  and  the  overall  relaxation  procedure  are  discussed. 


Inviscid  Flow  Representation 


To  accommodate  a  variable  mesh  into  the  difference  representation  of  the 
Cauchy  Integral,  Eq.  (17),  a  straightforward  extension  of  the  numerical  inte¬ 
gration  procedure  derived  by  Napolitano  et  al.  (Ref.  18)  has  been  introduced. 
The  slight  differences  that  are  encountered  are  due  solely  to  the  lack  of 
symmetry  of  the  variable  grid  with  respect  to  the  solution  point.  Thus,  in 
the  interaction  region  the  Cauchy  integral  is  approximated  by 


[p(Xi)-p,]  -■  r  Xj-Xj  1 


(63) 


1  "j  >'"[  1^1  ■<  Vi  ■’']>] 


where  the  function  h^  =  h  +  6  defines  the  displacement  surface,  2>,  the  sub¬ 
scripts  i  and  j  denote  streamwise  grid  point  indices,  and  I  is  the  index  of 
the  downstream  boundary  of  the  numerical  solution  domain. 


To  evaluate  the  contribution  to  the  local  pressure  from  displacement 
effects  upstream  of  the  strong  interaction  region;  i.e., 

.X,  h’^(0 

one  of  two  approaches  is  employed  depending  on  whether  the  flow  is  laminar 
or  turbulent.  For  the  laminar  case,  a  square-root  growth  in  displacement 
thickness  is  assumed  and  p^  is  evaluated  analytically.  For  turbulent  flow. 
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the  upstream  integral  is  evaluated  numerically  with  the  required  displace¬ 
ment  thickness  determined  from  a  direct  noninteracting  boundary  layer 
solution  for  the  given  airfoil  geometry. 


The  contribution  to  the  local  pressure  from  the  downstream  wake;  i.e.. 


p^(x)  =  - 


Li  x-C 


is  evaluated  using  the  asymptotic  expression 

S(x)  --  C,  +C2(x-Xor'/^ 


(66) 


for  the  far-wake  displacement  thickness,  which  holds  for  both  laminar  and 
turbulent  flow.  While  it  is  possible  to  relate  th^j  values  of  the  con¬ 
stants,  Cj^,  C2  and  Xq  to  the  airfoil  properties  (e.g.,  the  airfoil  drag 
determines  Cj^),  a  more  consistent  numerical  solution  is  achieved  by  eval¬ 
uating  these  constants  based  on  the  values  of  d  at  the  last  three  mesh 
points  of  the  strong- interact  ion  solution  region.  Thus,  the  constants 
in  Eq.  (66)  can  be  evaluated  using  the  following  difference  expressions 


Si-Si-,  +  2Si  Xi-2S;_,Xj., 


2(Sj-  Sj_,) 


“  8j  +  28j(xj~  Xq)  , 


Cg  ”  8j  28j  (xj“Xq) 


Viscous  Equation  Algorithm 

The  implicit  finite  difference  algorithm  for  Eqs.  (38)  and  (39)  differs 
only  slightly  from  that  presented  by  Werle  and  Verdon  (Ref.  4).  After 
introducing  the  computational  coordinates  S  =  S(£.)  and  N  =  N(C,n)>  Eq. 

(38)  can  be  written  as 

Vj=  Fj.,)+Qj 

(68) 


where  j  is  a  normal-direction  grid  point  index  and 

Pj  =  2(?jrT,  +  2RS  hj^/AS)/AN 

RS  an/as 


(69) 
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The  subscript  0  refers  to  the  adjacent  upstream  streamwise  station  and 

^  ^  d(lnS)  ^ 
d(lne)  ^ 
hjm  =  (hj.,  +  hj)/2 
~jm  =  (^j-i  +  n)/2 


Similarly,  the  momentum  Eq.  (39)  can  be  written  as 

Aj  Fj.,+  BjFj+  CjFj^,  +  DjVj  E,.  +  Eg. 


where 


Aj  =  ej^/(AN)^+  C^Vf./iAAU) 

Bj  =  -2RhjSCRF,./AS-(ejp+ej^)/(AN)^ 

Cj  =  ejp/{AN)^ -C^Vr. /(4AN) 

Dj  =  (CN/2-l)(Fr,^^-Fr._  V(2AN) 

E,.  =  Rhi(F.^-l) 

I]  J  rj 

Egj  =  -2RhjCRs[cL(2-CL)Foj  +  (l-CL)(2-CL)/(2F2j)]/AS 
+  CL(CL-l)FjFoj/2 
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liMji 


P  =  R/3 

~  ti'-i 

e  =  Kh 

Fr-  =  CL(2-CL)Foj  +  (l-CL)(2-CL)/(2F2j)  +  CL(CL-l)/(2Foj) 

(i-Cl)  Vgj 


0  implies  forward  linearization  and  convective  differencing 

Cj^  =  1  implies  backward  linearization  and  convective  differencing  ^ 

'2  implies  backward  iterative  linearization  and  convective 
differencing 


1  for  forward  flow 

Cr  = 

*'■  -1  for  reverse  flow 


The  coefficient  Cf^  refers  to  the  coupling  of  the  continuity  and  momentum 
equations  through  the  linearization  of  the  normal  convective  term. 


0  refers  to  the  first-order  linearization,  equations  coupled 


Cjj  =<1  refers  to  second-order  linearization,  equations  coupled 
12  refers  to  first-order  linearization,  equations  uncoupled 


*  These  values  were  incorrectly  given  in  reverse  order  in  Ref.  4 
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For  the  present  applications,  the  options  used  were  =  1,  Cj^  =  1,  and 
Cjj  =  1  except  in  reverse  flow  regions  where  Cl  and  were  set  equal  to 
0  and  -1,  respectively.  For  the  flat-plate  cases  where  a  singularity 
is  known  to  exist  at  the  trailing-edge  point,  the  difference  equations 
were  adjusted  for  a  first-order  coupled  linearization  (Cjj  =  0)  and  a  back¬ 
ward  convective  approximation  with  the  velocity  coefficient  linearized  about 
the  previous  iteration  level  (Cl  =2).  All  remaining  aspects  of  the  viscous 
solution  algorithm  are  identical  to  those  presented  in  Ref.  4. 


General  Solution  Algorithm 


The  computational  coordinates  are  defined  using  geometric  progres¬ 
sions  to  allow  a  packing  of  grid  points  in  regions  of  large  flow  gradients. 
In  the  normal  direction,  the  general  expression  used  for  setting  the  mesh 
spacing  and  hence  h  (n)  is 


1 


(74) 


where  K  is  the  mesh  spacing  constant,  is  the  value  of  the  turbulent 

Levy-Lees  variable,  n,  at  the  outer  edge  of  the  viscous  layer  and  AN  is  the 
normal  grid  spacing  in  the  computational  plane.  For  the  cases  considered 
here  forty  grid  points  have  been  used  in  the  normal  direction,  and  as  shown 
by  Carter,  Edwards,  and  Werle  (Ref.  19)  the  value  of  can  be  set  at  10 

for  virtually  all  flows  -  be  they  laminar,  transitional  or  turbulent. 

Selection  of  the  mesh  spacing  parameter,  K,  is  based  on  fluid  dynamic  require¬ 
ments.  For  laminar  flows,  K  is  prescribed  so  that  the  first  grid  point 
from  the  lower  boundary  is  at  the  same  physical  position  as  its  triple-deck 
solution  counterpart  (c.f.  Ref.  4).  It  was  found  that  this  could  be  accomp¬ 
lished  using  a  simple  uniform  spacing;  i.e.,  K  =  1.  For  turbulent  flow, 
mesh  strechlng  is  used  to  capture  the  anticipated  sub  or  inner  layer  of  the 
viscous  region.  It  was  found  that  using  a  value  of  K  =  1.05  provided  suf¬ 
ficient  resolution  of  the  inner  layer  for  the  turbulent  cases  considered  herein. 

Establishment  of  a  longitudinal  grid  spacing  involved  two  steps.  First,  the 
separation  coordinate  curve,  S  ,  is  positioned,  and  then  a  grid  distribution  is 
defined  using  a  double  geometric  progression  centered  at  the  trailing  edge. 

Since  only  first  derivatives  in  the  streamwise  direction  appear  in  the  govern¬ 
ing  equations,  it  is  only  necessary  to  provide  a  means  of  locating  the  indi¬ 
vidual  grid  points.  To  this  end  the  following  scheme  was  applied.  For 
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x<XjEi(l  =2,  ....  ijEi) 

Xj  =  X,  +  {Xt£-X,)(i- Ky'””)  /  (l-Ku 


(75a) 


3nci  for  ~  ^TEl*  ■  •  • » 


( i" i  j£ (  )  ( I~l) 

X:  =Xtp,  +{xtc,-Xte){i-K.  )/(i-K.  ) 


(75b) 


where  the  subscripts  TE  and  TEl  refer  to  grid  point  locations  at  the 
trailing  edge  and  immediately  downstream  of  the  trailing  edge,  respectively. 
The  parameters  used  to  adjust  the  upstream  grid  spacing  are  the  spacing 
parameter,  K^,  the  number  of  points,  i-pg,  and  the  overall  length,  -  Xj^ , 

For  the  downstream  region,  the  only  parameters  available  are  the  spacing 
constant,  K^,  and  the  total  number  of  points,  I,  since  here  the  overall 
length  is  established  from  the  value  of  -  x^^  which  is  set  by  the  up¬ 
stream  expression.  For  the  current  calculations  x^  was  generally  set  equal 

to  0.4  with  41  grid  points  ahead  of  the  trailing  edge.  Downstream  41  points 

were  again  used  over  a  distance  of  approximately  0.6.  For  laminar  flows, 
grid  stretching  was  not  applied  in  the  streamwise  direction  while  for  several 
of  the  turbulent  flow  cases,  Ky  and  were  set  equal  to  0.95  and  1.05, 
respectively.  For  blunt  trailing-edge  geometries,  the  shear-layer  coor¬ 
dinates  (Fig.  1)  are  employed.  The  mesh  stretching  is  then  centered  around 
the  shear-layer  coordinate  trailing-edge  point  which  is  located  at  the  inter¬ 
section  of  the  shear-layer  coordinate  curve  and  the  normal  coordinate  line 
which  passes  through  the  trailing  edge  point. 

As  depicted  in  the  sketch  below,  the  shear-layer  coordinate  curve  was  taken 
as  a  three  piece  fit  with  the  airfoil  surface  forming  the  initial  segment 
(xiX2),  followed  by  a  cubic  section  for  X2<x<x^  which  is  tangent  to  the  air¬ 
foil  surface  at  the  point  X2  and  to  the  wake  centerline  at  the  point  and 

with  the  wake  centerline  forming  the  third  segment  for  x>x^ .  The  trailing- 
edge  region  cubic  is  established  such  that  the  trailing-edge  location  on  the 
s-coordlnate  curve  is  positioned  at  a  prescribed  ratio,  r,  of  the  distance 
from  X2  to  That  is,  x^  is  set  according  to  the  relation 


r  =  (X3- X2)/(X4  -  Xg) 


(76) 


and  the  value  of  x^  is  computed  for  a  given  cubic  using  the  relation 
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X3-X2 


■  T  'y(’‘TE~’^2^^'‘’24h(K2)r^(i-r)  [i  +  r2(2r-3)]  ] 


(77) 


Here  r  was  set  equal  to  0.3  for  all  cases  studied. 


The  s-coordinate  curve  is  then  positioned  according  to  the  relation 


(78) 


and  the  coordinate  surface  inclination  angle  is  given  by 

6h(x2)  /  X-X2  \  /  X-X2 


(79) 


Finally,  there  is  a  need  to  establish  the  airfoil  height,  h(xjj),  at  a  coor¬ 
dinate  location  s^,;  i.e.. 


h(xh)  =  [b^-  2a(xb-XTE)- c  (Xb-x.^)^]'^^-b 


(80a) 
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where 


"  ^TE 


(3  -7i  tan^  +  [(d  a  tor?  -  U  +  c  {a'^-b^tan^^;^.)]' 


(80b) 


(80c) 


1  +  C 


a  =eT£(ho/eTEJ  7/(y-2) 

b  =  ho/(x-2) 

C  -  (h/6^/  y/{y-2) 
d  =  xc~xte+  {h(xc)  +  b)  ton^j. 
y  =  e-TE^Q^af/ho 


To  apply  the  present  iterative  solution  procedure,  it  is  necessary  to 
establish  initial  estimates  to  the  displacement  surface  and  the  boundary  layer 
profile  at  the  upstream  boundary,  ^  of  the  computational  domain.  This 

information  was  determined  using  the  viscous  solution  algorithm  in  the  direct 
mode  which  provides  the  required  velocity  profile  and  displacement  thickness 
distribution  for  a  specified  pressure  distribution.  This  imposed  pressure 
distribution  was  calculated  using  the  inviscid  pressure  Eq.  (17),  with  the 
shear-layer  coordinate  curve  taken  as  the  first  estimate  to  the  displacement 
surface.  With  the  resulting  pressure  distribution,  a  marching  viscous  solu¬ 
tion  proceeds  from  the  leading  edge  of  the  airfoil  to  the  upstream  boundary, 

X  =  gf  the  strong- interaction  region  to  establish  the  initial  velocity 
profile  and  the  displacement  contribution  to  the  upstream  pressure  integral,  Eq. 
(64).  For  a  flat-plate  airfoil  this  initial  marching  solution  is  unnecessary 
since  a  self-similar  boundary  layer  solution  can  be  determined. 

Subsequent  interative  passes  through  the  inviscid  and  viscous  solutions 
proceeded  with  a  relaxation  factor,  w,  generally  held  at  a  value  near  unity, 
to  achieve  a  converged  solution  in  which  the  displacement  surface  ceased  to 
vary  in  the  4th  significant  figure.  For  some  of  the  larger  separation  bubble 
cases,  this  level  of  convergence  required  approximately  50  Iterations,  with 
a  relaxation  factor  of  1/3. 
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APPLICATION  TO  TRAILING-EDGE  FLOWS  -  RESULTS  AND  DISCUSSION 


The  analytical  and  numerical  approach  outlined  above  has  been  assessed 
through  a  series  of  model  problem  studies.  The  first  series  focused  on 
laminar  flat-plate  traillng-edge  flow  to  provide  a  fundamental  assessment 
of  the  present  approach  through  detailed  comparisons  with  triple-deck  solu¬ 
tions.  The  second  test  series  considered  the  extension  of  the  interacting 
boundary  layer  (IBL)  approach  to  laminar  blunt-based  configurations  where 
again  triple  deck  solutions  are  available  for  comparison.  The  third  study 
involved  the  introduction  of  turbulence  effects  with  initial  focus  on  flat- 
plate  trailing-edge  flow  where  limited  data  is  available  for  assessment.  The 
final  test  series  considered  application  of  the  interacting  boundary  layer 
approach  to  turbulent  blunt-based  f lows-lncluding  separation  phenomena. 


Laminar  Flow 


Flat-Plate  Triple-Deck  Solutions 

As  discussed  in  Ref.  4,  the  basic  numerical  algorithm  employed  in  the 
interacting  boundary  layer  analysis  is  a  direct  extension  of  the  numerical 
method  developed  in  Ref.  4  for  the  solution  of  the  triple-deck,  asymptotic- 
limit  (Re-«“),  form  of  the  governing  equations.  It  is  thus  imperative  that 
this  triple-deck  algorithm  first  be  assessed  by  comparison  of  results  with 
those  of  other  investigators  (e.g..  Refs.  20  and  21).  The  problem  of  calcu¬ 
lating  the  laminar  solution  for  the  trailing-edge  region  of  a  flat  plate 
imposes  a  severe  test  to  a  numerical  solution  algorithm  because  of  the 
occurrence  of  a  fundamental  singularity  in  the  pressure  gradient  at  the 
trailing  edge.  Current  results,  obtained  with  the  basic  algorithm  (desig¬ 
nated  by  the  solid  curve  in  Fig.  2)  compare  very  favorably  with  those  of 
Refs,  20  and  21  up  to  the  trailing  edge  but  significant  differences  occur 
along  the  wake  (c.f..  Fig.  2).  Note  that  the  pressure,  P,  and  distance,  X, 
(Fig.  2)  are  dimensionless  variables  which  have  been  scaled  according  to 
triple-deck  scaling  laws.  There  are  numerous  possible  sources  for  this 
discrepancy,  ranging  from  truncation  errors  to  an  inaccurate  representation 
of  the  trailing-edge  singularity  (see  Ref.  20  for  a  detailed  discussion  of 
this).  Of  these,  it  was  found  that  the  most  influential  was  the  manner  in 
which  the  viscous  pressure  distribution  was  determined.  If  is  obtained 

by  a  forward-marching  integration  across  the  mesh  element  just  aft  of  the 
trailing  edge,  the  singularity  at  the  trailing  edge  will  not  be  properly 
accounted  for.  To  avoid  this  error,  solutions  have  been  obtained  with  a 
forward-marching  Integration  up  to  the  trailing  edge  and  a  backward- 
marching  Integration  (from  the  downstream  asymptotic  solution)  to  the  mesh 
element  just  aft  of  the  trailing  edge.  Wake  pressure  distributions 
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determined  with  this  new  integration  scheme  (Fig,  2(a))  are  in  better  agreement 
with  those  of  previous  investigators.  The  remaining  differences  are  attributed 
to  truncation  error,  as  demonstrated  in  Fig,  2(b),  which  shows  the  approach  of 
the  current  trailing-edge  region  pressure  distributions  toward  the  analytic 
predictions  of  Ref.  20  as  the  longitudinal  grid  spacing  is  reduced  from  AX  = 

0.2  to  0.05.  It  is  important  to  note  that  in  both  Refs.  20  and  21,  special 
analytical  adjustments  were  applied  at  the  trailing  edge. 

Flat-Plate,  Finite-Reynolds-Number  Solutions 

Interacting  boundary  layer  solutions  have  been  determined  to  verify 
their  systematic  approach  to  the  triple-deck  limit  solution  as  Re-+“.  For 
this  purpose,  two  separate  studies  have  been  conducted.  In  the  first,  the 
forward/reverse  viscous  pressure  integration  discussed  relative  to  Fig.  2 
was  employed  and  in  the  second,  the  forward -marching  integration  was  used. 

The  resulting  comparisons  of  the  IBL  solutions  with  their  respective  triple¬ 
deck  solutions  are  shown  in  Figs.  3(a)  and  3(b)  for  a  longitudinal  step  size 
of  AX  =  0.2.  The  IBL  solutions  were  obtained  with  the  physical  mesh  spacing 
adjusted  so  that  AX  =  0.2  and  the  variable  normal  grid  spacing  set  such  that 
the  spacing  in  the  inner  region  of  the  viscous  layer  corresponded  to  AY  =  0.4, 
where  Y  is  the  triple-deck  normal  coordinate.  As  shown  in  Figs.  3(a)  and 
3(b),  the  resulting  IBL  solutions  shov;  a  systematic  approach  to  their  respec¬ 
tive  triple-deck  counterparts.  Iifhile  in  an  overall  sense,  these  results 
verify  the  IBL  algorithm,  they  also  indicate  that  the  trailing-edge  singular¬ 
ity  is  causing  complications  that  must  be  addressed.  The  Davis  (Ref.  22) 
algorithm  apparently  avoids  these  complications  and  consideration  should  be 
given  to  adopting  several  features  of  his  approach  in  our  future  studies. 

One  aspect  of  the  interacting  boundary  layer  solution  which  was  found 
to  be  of  paramount  importance  is  the  contribution  to  the  pressure  in  the 
strong-interaction  region  from  displacement  effects  upstream  and  downstream 
of  this  region  (c.f.,  Eqs.  (64)  and  (65)).  To  illustrate  this  point,  solu¬ 
tions  have  first  been  determined  for  a  flat  plate  of  infinite  extent  where 
formal  weak- interact ion  theory  predicts  that  the  pressure  coefficient, 

Cp  =  P-P„,  is  zero  over  the  length  of  the  plate.  Solutions  obtained  for 
Re  =  105  with  the  computational  region  extending  from  s  =  x  =  0.4  to  s  =  x  = 

1.6  are  shown  in  Fig.  3(c).  A  comparison  of  the  solutions  with  and  without 
the  far-field  contributions,  p.^  and  p^,  to  the  local  pressure  clearly  shows 
the  significant  impact  of  these  terms  at  this  Reynolds  number.  Figure  3(c) 
also  compares  solutions  for  a  finite  flat  plate  obtained  with  and  without 
the  far-field  contributions  to  the  local  pressure.  The  results  indicate 
rather  large  differences  especially  near  the  computational  boundaries. 

Clearly,  far-field,  weak-interaction  phenomena  can  be  expected  to  have  a 
significant  impact  on  solutions  determined  in  local  strong-interaction  regions. 
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Solutions  to  the  interacting  boundary  layer  equations  for  a  finite  flat 
plate  at  a  Reynolds  number  of  10^  have  also  been  presented  by  Davis  (Ref.  22) 
and  Veldman  (Refs.  23-25).  Their  results  are  compared  with  the  present  solu¬ 
tion  in  Figs.  4(a)  through  4(c).  Veldman’s  solutions,  based  on  a  composite 
boundary-layer/triple-deck  model  and  on  the  interacting  boundary  layer  model 
were  obtained  without  accounting  for  far-field  displacement  effects  on  the 
upstream  flow  as  discussed  above.  Thus,  it  is  not  surprising  that  consider¬ 
able  differences  exist  between  predicted  upstream  pressure  levels  (Fig.  4(a)). 
However,  present  results  compare  very  favorably  with  the  recent  results  of 
Davis  (Ref.  22)  which  were  obtained  with  global  far-field  effects  accounted 
for.  Only  minor  differences  (Fig.  4(a))  occur  just  aft  of  the  trailing 
edge  where  very  steep  gradients  are  encountered.  Based  on  the  mesh  size 
study  described  above,  one  would  expect  these  differences  to  vanish  with 
decreasing  mesh  spacing. 

Results  for  skin  friction,  wake  centerline  velocity,  and  displacement 
thickness  distributions  are  presented  in  Figs.  4(b)  and  4(c).  The  various 
interaction  solutions  are  seen  to  be  in  excellent  agreement.  Also  shown 
are  the  "noninteracting"  results  obtained  from  the  boundary  layer/wake 
equations  with  the  imposed  Inviscld  pressure  gradient  taken  to  be  zero. 

The  noninteracting  displacement  thickness  distribution,  Fig.  4(c),  has  a 
square-root  growth  up  to  the  trailing  edge  followed  by  a  rapid  (singular) 
drop-off  in  the  wake  just  aft  of  the  trailing  edge,  which  is  due  to  the 
change  in  surface  boundary  condition  at  the  trailing  edge.  The  singularity 
at  the  trailing  edge  is  relieved  when  strong-interaction  effects  are 
accounted  for  (IBL  analysis).  It  can  be  seen  that  the  strong  interaction 
at  the  trailing  edge  is  felt  upstream  up  to  approximately  x  =  0.75  and  that 
it  significantly  alters  the  displacement  thickness  gradient  in  the  near  wake 
(1.0^  X  <1.25).  Seen  in  this  light,  the  strong-interaction  region  appears 
to  be  confined  to  the  region  0.75<x<1.25  where  it  provides  a  local  correction 
to  the  global  near-wake  neck-down  effect.  This  strong  interaction  causes 
an  acceleration  and  thinning  of  the  boundary  layer  near  the  trailing  edge 
(Fig.  4(c))  with  a  subsequent  increase  in  surface  shear  stress  (Fig.  4(b)). 

Blunt  Traillng-Edge  Solutions 

The  present  numerical  scheme  has  been  applied  to  blunt  trailing-edge  flows 
typical  of  those  encountered  in  turbomachinery  applications.  A  shear-layer  coordi¬ 
nate  system  has  been  established  and  calculations  have  been  performed  for  laminar 
flows  past  thickened  flat-plate  airfoils  (oy,  =  0)  with  elliptic-section  trailing 
edges  (see  Fig.  5(ci)).  Converged  solutions  have  been  obtained  for  attached  and 
separated  flows  at  Reynolds  numbers  (based  on  plate  length)  of  10^  and  10^  which 
will  bj  shown  to  approach  their  triple-deck  limits  with  increasing  Re  for  attached 
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mildly  separated  flows.  Sample  results  for  a  free-stream  Reynolds  number  of 
105  and  Mach  number  of  0.1  with  an  elliptic  section  length,  ,  of  0.2  and 
heights,  h^,  of  0.07  and  0.14  are  shovm  in  Figs.  5(a)  through  5(f).  These 
parameters  give  an  inner  deck  scaling  of  L  =  4,  Hq  =  4  and  8,  and  y  =  2, 
and  thus  correspond  to  two  of  the  triple-deck  solutions  reported  in  Ref.  4. 
The  interacting  boundary  layer  solutions  were  obtained  for  a  strong-inter¬ 
action  calculation  region  extending  from  40  percent  chord  to  160  percent 
chord.  Classical  Blasius  boundary  layers  and  Goldstein  wakes  were  assumed 
to  exist  upstream  and  downstream,  respectively,  of  the  Interaction  region. 

The  shear-layer  coordinate  curve  was  taken  as  a  simple  cubic  (see  Fig.  5(b)), 
tangent  to  the  airfoil  at  the  trailing-edge  juncture  point  and  to  the  wake 
centerline  at  140  percent  chord. 

In  the  resulting  solutions  (Figs.  5(b)  through  5(f))  the  displacement 
surfaces  show  steady  growth  up  to  near  the  juncture  point  followed  by  a 
rapid  smooth  decay  in  the  near  wake  toward  constant  far-wake  values 
(Fig.  5(b)).  The  surface  skin  friction  and  wake  centerline  velocity  dis¬ 
tributions  are  shown  in  Fig.  5(c)  where  the  thick  blade  case,  =  8,  is 
seen  to  experience  considerable  separation  off  the  blunt  base.  The  oscilla¬ 
tions  in  the  skin  friction  at  the  upstream  boundary  of  the  interaction 
region  are  due  to  the  rapid  adjustment  of  the  assumed  initial  self-similar 
velocity  profile  to  a  compatible  nonsimilar  state.  The  thick  plate  case 
experiences  a  skin  friction  growth  Immediately  downstream  of  the  initial 
station  which  is  due  to  boundary  layer  thinning  induced  by  a  flow  overspeed 
as  it  accelerates  around  the  elliptic  trailing  edge. 

The  displacement  thickness  distribution  over  the  trailing  edge  region 
(c.f..  Fig.  5(d))  has  a  singularity  at  the  trailing  edge  point.  If  the 
shear-layer  coordinate  curve  (c.f..  Fig.  5(b))  were  taken  to  coincide 
with  the  airfoil  surface;  i.e.,  t  =  0,  then  the  classical  first-order 
Invlscid  flow  would  be  singular  at  the  trailing  edge.  The  current  approach 
eliminates  this  singularity  by  defining  a  smooth  coordinate  surface  in 
such  a  fashion  that  the  blunt  trailing-edge  singularity  is  removed  from  the 
first-order  Invlscid  solution.  It  should  be  noted  that  the  shear-layer 
coordinate  concept  has  been  recently  applied  by  Carter  et  al.,  (Ref.  26) 
for  solution  of  the  transonic  shock-induced  separation  problem.  This 
application  represents  one  of  the  major  contributions  of  the  current  effect. 

The  pressure  distributions  for  the  blunt  trailing-edge  geometries  are 
shown  in  Figs.  5(e)  and  5(f).  For  these  cases,  the  pressure  levels  from 
the  successive  invlscid  and  viscous  calculations  converge  relatively  slowly 
in  the  near-wake  region.  This  is  illustrated  in  Fig.  5(f)  where  the  respec¬ 
tive  Inviscld  and  viscous  pressure  levels  are  shown  after  50  iterations. 
However,  it  was  found  that  the  mean  pressure  distribution  (as  shown  in 
Fig.  5(e))  was  well  converged  after  only  20  iterations  and  that  both  the 


32 


R81-914986-5 


Invlscld  and  viscous  pressures  headed  toward  this  level.  This  result 
Indicates  that  the  Iteration  scheme  might  be  Improved  through  use  of  an 
averaging  technique. 

The  pressure  distributions  depicted  in  Fig.  5(e)  show  Initial  pressure 
levels  significantly  different  from  their  flat-plate  counterpart 
(see  Fig.  4(a)).  This  is  due  both  to  the  blunt  geometry  and  global  weak- 
interaction  effects.  As  discussed  above,  the  initial  pressure  level  is 
established  by  adding  the  global  weak-interaction  contribution  (i.e.,  the 
far-fleld  contribution  to  the  pressure  in  the  strong-interaction  region) 
to  the  geometric-induced  pressure  level.  The  pressure  drop  downstream  of 
the  initial  station  is  caused  by  the  flow  acceleration  around  the  blunt 
trailing  edge.  There  is  a  rapid  pressure  rise  across  the  trailing  edge  to 
a  near-wake  peak  value  followed  by  a  weak-interaction  pressure  decay  in  the 
far  wake.  The  appearance  of  a  "plateau"  like  region  for  the  separated 
case,  =  8,  leads  to  a  pressure  signature  in  the  near  wake  region  quite 
similar  to  that  observed  in  experimental  studies  of  bluff-based  bodies. 

It  is  surprising  to  find  a  very  fine  scale  pressure-gradient  reversal  at 
trailing  edge  for  the  separated  case  and  it  has  not  been  established 
whether  this  is  realistic  or  merely  a  feature  of  the  present  interaction 
model. 

A  final  assessment  of  the  interacting  boundary  layer  code  for  blunt 
trailing  edges  is  obtained  through  comparison  with  the  triple  deck  solutions 
of  Ref.  4.  The  only  available  basis  for  assessing  the  current  blunt  airfoil 
solutions  is  to  ascertain  whether  they  approach  corresponding  triple-deck 
solutions  as  the  Reynolds  number  is  increased.  To  perform  such  a  study  it 
is  necessary  to  adjust  the  scale  of  the  physical  geometry  and  finite  differ¬ 
ence  mesh  such  that  they  honor  the  triple-deck  scaling  laws.  Figures  6(a) 
through  6(c)  present  the  results  of  such  a  study  for  two  thicknesses,  = 

4  and  H  =  8,  as  Re-«“.  The  results  have  been  presented  in  triple-deck 
variables  with  the  asymptotic  solutions  (Re-«°)  taken  directly  from  Ref.  4. 

In  general,  the  results  show  a  relatively  insensitive  dependence  on  Reynolds 
number  Indicating  that  the  asymptotic  solution,  Re-»-“>,  can  be  used  to  obtain 
good  qualitative  and  quantitative  predictions  for  a  Reynolds  number  as  low 
as  10^.  This  result  is  in  sharp  contrast  to  the  supersonic  separated  flow 
results  obtained  by  Burggraf  et  al.,  (Ref.  27)  where  the  density  variations 
across  the  inner  layer  severely  limited  the  range  of  applicability  of  triple- 
deck  theory. 

Considering  the  specific  results  of  this  study,  attention  is  first 
directed  to  the  unseparated  case,  =  4.  The  IBL  results  show  a  slow  but 
systematic  approach  to  the  triple-deck  solution  over  the  entire  solution 
domain  (Fig.  6(a)).  The  accompanying  shear  stress  and  wake  centerline 
velocity  distributions  (Fig.  6(b))  behave  similarly  except  near  the  Initial 
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station  (X  =  -6)  of  the  computational  region.  Such  behavior  is  believed  to 
be  traceable  to  the  initialization  procedure.  For  the  IBL  results,  the  mesh 
spacing  was  determined  such  that  in  triple-deck  variables  it  was  fixed  at 
AX  =  0.1  as  Re-x®.  For  the  fixed  number  of  longitudinal  mesh  points  used 
in  the  current  study  (121),  this  corresponded  to  an  upstream  boundary  at 
X  =  -6.  Apparently  this  location  is  too  close  to  the  trailing  edge  to 
allow  the  use  of  a  self-similar  velocity  profile  to  start  the  calculation. 
Fortunately,  initial  errors  appear  to  damp  out  rather  quickly. 

The  situation  is  slightly  changed  for  the  thicker  airfoil,  =  8. 

Again  there  is  an  anamoly  in  shear  stress  near  the  initial  station,  X  =  -6, 

(Fig.  6(c))  but  this  effect  quickly  disappears  as  X  increases.  A 

notable  difference  occurs  in  the  wake  region  where  the  pressure  (Fig.  6(a)) 

and  wake-centerline  velocity  (Fig.  6(c))  distributions  show  a  small  but  distinct 

sensitivity  to  Re  in  the  near  wake.  Note  that,  as  in  the  Hq  =  4  case,  the  starting 

procedure  apparently  causes  a  nonuniformity  in  the  behavior  of  pressure  and  skin 

friction  with  increasing  Reynolds  number  over  the  interval  -4<X<-1  (as  shown  in 

Figs.  6(a)  and  6(c)). 


Turbulent  Flow 

Application  of  the  current  analysis  to  the  case  of  turbulent  trailing- 
edge  flows  was  accomplished  using  the  eddy  viscosity  formulation  outlined  in 
Eqs,  (28)  through  (30).  This  eddy  viscosity  model  has  been  expressed  in 
terms  of  the  "turbulent"  Levy-Lees  type  variables,  Eqs.  (51)  through  (53), 
in  which  the  longitudinal  independent  variable  is  defined  in  such  a  manner 
as  to  minimize  the  effects  of  turbulence  on  the  growth  of  the  viscous  region 
in  the  computational  domain.  In  addition,  the  Levy-Lees  normal  coordinate 
has  been  scaled  to  provide  more  mesh  points  near  the  wall  and  wake  center- 
line  so  that  a  better  resolution  of  the  turbulent  inner  layer  and  also  of 
thin  separation  bubbles  can  be  obtained. 

For  the  turbulent  cases,  the  pressure  contribution  from  the  region 
ahead  of  the  strong-interaction  zone  (Eq.  (64))  cannot  be  evaluated  analyt¬ 
ically.  Thus  the  upstream  displacement  thickness  distribution  was  taken 
directly  from  classical  non-interacting  boundary  layer  calculations  with 
instantaneous  transition  Imposed  immediately  aft  of  the  leading  edge.  The 
appropriate  far-wake  displacement  thickness  distribution  is  of  the  Goldstein 
square-root  decay  type  and  thus  the  procedure  employed  for  the  laminar  case 
carried  over  directly. 

Flat-Plate  Trailing  Edge 

After  the  turbulence  model  and  the  variable  mesh  spacing  capability 
were  incorporated  into  the  finite  Re  interaction  algorithm,  a  viscous/ 
Inviscld  interaction  solution  was  determined  for  subsonic  turbulent  flow 
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at  Moo  ■  0.0116  and  Re  ■  6.5  x  10^  past  a  flat-plate  trailing  edge.  The 
vlscous/inviscld  interaction  pressure  distribution  for  this  case  is  shown 
in  Fig.  7(a)  along  with  a  prediction  based  on  a  second-order,  noninteracting, 
boundary-layer  solution.  The  noninteracting  solution  was  determined  over 
the  length  of  the  plate  and  into  the  wake  region  to  obtain  a  first-order 
displacement  thickness  distribution.  In  the  spirit  of  higher-order  boundary 
layer  corrections,  the  inviscid  pressure  distribution  was  then  calculated 
based  on  this  displacement  thickness  distribution.  The  result  shown  in 
Fig.  7(a)  reveals  an  expected  singular  behavior  at  the  trailing  edge.  The 
interaction  solution  indicates  that  the  displacement  interaction  effect  is 
highly  confined  to  near  the  trailing  edge  point,  as  anticipated  by  Melnik 
(Ref.  1),  and  serves  principally  to  remove  the  pressure  discontinuity  at 
the  trailing  edge.  The  conclusion  drawn  here  is  that  while  strong  inter¬ 
action  effects  are  important,  they  do  not  dominate  weak,  global  interaction 
effects  except  in  the  immediate  vicinity  of  the  trailing  edge. 

Present  results  for  wake  centerline  velocity  are  compared  with  data 
taken  by  Chevray  and  Kovasznay  (Ref.  28)  and  results  of  the  noninteracting 
solution  of  Cebeci  et.al.  (Ref.  15)  in  Fig.  7  (b).  It  can  be  seen  that  the 
agreement  between  the  noninteracting  and  interacting  solutions  and  the  data 
is  very  good.  Important  features  are  that  the  interacting  solution  gives  a 
slight  Improvement  over  the  noninteracting  results  in  the  near  wake,  and 
that  both  solutions  underpredict  the  wake-centerline  velocity  in  the  far 
wake.  Such  differences  are  well  within  the  uncertainty  surrounding  the 
turbulence  model  for  trailing-edge  flow. 

Blunt  Trailing  Edge 

The  final  cases  presented  here  are  those  for  turbulent  flow  off  the 
trailing  edge  of  a  thickened  flat  plate  airfoil.  The  geometry,  Reynolds 
number,  and  Mach  number  are  identical  to  those  studied  for  laminar  flow 
(c.f..  Figs.  5(a)  and  6).  The  resulting  displacement  body  and  surface 
pressure  distributions  for  =  A  are  compared  to  their  laminar  counter¬ 
parts  in  Figs.  8(a)  through  8(c).  The  turbulent  results  show  steeper 
gradients  in  displacement  thickness  and  pressure  across  the  trailing-edge 
region  than  those  encountered  in  the  laminar  case.  This  results  in  a 
slightly  higher  peak  pressure  in  the  near  wake.  The  far-wake  displacement 
thickness  is  larger  for  turbulent  than  for  laminar  flow,  which  is  apparently 
due  to  the  increased  drag  induced  on  the  airfoil  in  the  turbulent  case. 

It  is  interesting  to  note  the  degree  to  which  the  local  interaction  effect 
Influences  the  pressure  distribution.  The  rapid  decrease  in  displacement 
thickness  in  the  near  wake  so  dominates  the  problem  that  significant  dif¬ 
ferences  between  laminar  and  turbulent  pressure  distributions  can  only  be 
detected  over  a  relatively  small  region  of  the  near  wake. 
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Results  indicating  the  effect  of  airfoil  thickness  on  turbulent 
trailing-edge  flows  are  shown  in  Figs.  9(a)  through  9(c).  The  strong  inter¬ 
action  solutions  were  initiated  at  a  40  percent  chordwise  station  with  the 
initial  boundary  layer  profiles  determined  from  noninteracting  solutions 
upstream  of  this  region.  Minor  nonuniformities  appear  in  the  skin  friction 
parameter  and  pressure  distributions  of  Figs.  9(b)  and  9(c)  near  the  initial 
station  because  trailing-edge  interaction  effects  apparently  do  not  fully 
attenuate  at  this  upstream  location.  The  far-wake  displacement  thickness 
was  again  represented  asymptotically  as  a  square-root  decay  to  a  constant 
level.  The  constants  of  this  asymptotic  representation  were  determined 
by  the  strong-interaction  solution.  As  shown  in  Fig.  9(a),  the  shear-layer 
coordinate  curve  was  set  tangent  to  the  airfoil  surface  at  the  trailing-edge 
juncture  point  and  to  the  wake  centerline  at  125  percent  chord. 

The  displacement  surface  predictions  are  shown  in  Fig.  9(a)  where  a 
steady  growth  is  indicated  up  to  the  trailing-edge  juncture  point,  followed 
by  a  very  rapid  thinning  into  the  near  wake,  and  then  a  gradual  thinning 
toward  constant  levels  in  the  far  wake.  Note  that  the  distorted  scale  of 
Fig.  9(a)  tends  to  over-dramatize  this  effect  and  that  true  flow  angles  are 
much  less  severe  than  those  shown.  It  is  surprising  to  find  that  the  far- 
wake  displacement  body  approaches  the  same  constant  level  for  the  two 
cases  even  though  these  displacement  bodies  differ  by  as  much  as  70  percent 
at  the  strong- interaction  initiation  point.  This  can  be  explained  to  some 
extent  by  a  consideration  of  the  asymptotic  wake  structure.  The  wake  dis¬ 
placement  thickness  is  a  direct  function  of  airfoil  drag.  Since  tlie  airfoils 
are  of  similar  geometry,  one  would  expect  and  Fig.  9(b)  verifies,  that  the 
skin  friction  drag  would  he  nearly  identical  for  the  two  cases.  The  only 
difference  then  must  be  traceable  to  the  Immediate  trailing-edge  region, 
which  is  physically  of  very  small  extent. 

The  skin  friction  parameter  and  wake  centerline  velocity  distributions 
of  Fig.  9(b)  indicate  very  small  differences  over  the  entire  trailing-edge 
region  with  only  a  very  mild  separation  occurring  for  the  thick  airfoil 
(Hq  =  8).  However,  the  pressure  distributions  of  Fig.  9(c)  do  show  inter¬ 
esting  and  significant  differences.  There  are  three  major  contributions 
to  this  pressure  difference;  first,  that  due  to  geometric  differences  which 
sets  overall  lower  surface  and  higher  wake  pressures  for  the  thick  airfoil; 
second,  the  global  weak  interaction  contribution  from  the  finite-thickness, 
far-wake  displacement  body  which  mitigates  the  geometric  differences;  and 
finally,  the  near-wake  strong- interaction  effects  whicli  serve  to  bridge  the 
differences  between  the  respective  surface  and  wake  pressure  distributions. 
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CONCLUSIONS  AND  RECOMMENDATIONS 


A  general  approach  based  on  viscous/lnviscld  interaction  theory  has  been 
outlined  for  determining  high  Reynolds  number  flow  past  a  cascade  of  airfoils. 
The  overall  flow  field  is  broken  down  into  "inner"  dissipative  fields  consist¬ 
ing  of  attached  boundary  layers,  wakes,  and  local  strong-interaction  (possibly 
separated)  flow  regions  and  an  "outer"  inviscid  field.  The  flow  in  the  inner 
or  viscous  regions  is  assumed  to  be  governed  by  Prandtl's  boundary  layer  equa¬ 
tions  and  the  outer  flow  is  assumed  to  be  isentropic  and  irrotational  (i.e., 
the  Interacting  Boundary  Layer  Model).  Solutions  to  the  viscous  and  inviscid 
equations  are  matched  to  provide  a  solution  for  the  entire  flow  field.  The 
matching  procedure  depends  on  the  nature  (i.e.,  weak  or  strong)  of  the  inter¬ 
action  between  the  inviscid  and  viscous  flows. 

In  the  present  study  emphasis  has  been  placed  on  the  development  of  a 
solution  procedure  for  local  strong-interaction  regions  and,  in  particular, 
the  strong-interaction  region  caused  by  the  change  in  boundary  condition  exper¬ 
ienced  by  the  fluid  as  it  leaves  the  surface  at  the  trailing  edge  of  an  air¬ 
foil.  Subsonic  laminar  and  turbulent  flows  past  sharp  and  blunt  trailing  edges 
have  been  considered.  A  shear-layer  coordinate,  which  departs  from  the  airfoil 
surface  upstream  of  the  trailing  edge  and  "reattaches"  on  the  wake  centerline, 
has  been  Introduced  for  blunt-edge  applications.  After  making  the  usual  thin 
viscous  layer  assumptions  relative  to  the  shear-layer  coordinate  system  and 
applying  Prandtl's  transposition  theorem,  the  classical  form  of  the  boundary 
layer  equations  is  retained.  These  equations  are  then  recast  in  terms  of 
"turbulent"  Levy-Lees  variables  to  facilitate  numerical  resolution.  For  strong 
viscous/inviscid  interactions  an  iterative  solution  procedure  is  required.  In 
the  present  effort  successive  inviscid  and  viscous  solutions  are  determined 
directly  and  inversely,  respectively;  i.e.,  the  displacement  thickness  is  pre¬ 
scribed  and  the  two  governing  equation  sets  are  solved  to  determine  the  inviscid 
and  viscous  pressure  fields.  Updated  values  of  displacement  thickness  are  ob¬ 
tained  using  Carter's  (Ref.  6)  relaxation  procedure  until  a  converged  solution 
for  the  overall  flow  field  is  achieved. 

Numerical  results  have  been  presented  here  for  symmetric  adiabatic  flows 
past  thln-alrfoll  trailing  edges.  Laminar  flat-plate  results  are  in  very  good 
agreement  with  those  of  other  Investigators  (Refs.  22  and  23).  In  addition,  the 
present  finite  Reynolds  number  solutions  approach  asymptotic  (Re  “)  triple¬ 
deck  solutions  (Refs.  20  and  21)  with  increasing  Reynolds  number.  Separated 
laminar  flow  solutions  have  been  determined  for  thickened  flat-plate  airfoils 
with  elliptic-section  trailing  edges.  These  blunt-edged  results  also  approach 
triple-deck  solutions  (Refs.  A  and  5)  with  increasing  Reynolds  number.  Turbu¬ 
lent  flows  past  flat-plate  and  thickened  flat-plate  airfoils  have  also  been 
calculated.  For  these  cases,  no  other  interaction  solutions  are  available  for 
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assessing  the  current  analysis.  However,  present  results  for  wake-centerline 
velocities  behind  a  thin  flat-plate  airfoil  have  been  found  to  be  in  good 
agreement  with  the  experimental  data  reported  by  Chevray  and  Kovasznay 
(Ref.  28). 

In  the  present  effort,  consideration  has  been  restricted  to  symmetric  flows 
and  hence,  a  primary  goal  in  future  work  will  be  the  extension  of  the  present 
concepts  and  solution  procedures  to  asymmetric  traiiing-edge  flows  which  are 
representative  of  those  encountered  in  turbomachinery  applications.  In  addi¬ 
tion,  there  is  a  need  for  assessing  the  current  analysis  through  comparison 
with  experimental  measurements  of  pressure  and  velocity  distributions  in  the 
turbulent  near-wake  region  of  a  symmetric  airfoil.  Finally,  the  existing 
numerical  solution  procedure  should  be  further  developed  to  improve  the  resol a- 
tlon  of  flows  with  relatively  large  separation  regions,  and  to  improve  the 
convergence  rate  of  Che  viscous/inviscid  iteration  procedure.  One  possibility 
here  is  to  incorporate  the  quasi-simultaneous  iteration  procedure  developed 
by  Veldman  (Ref.  8)  and  Davis  (Ref.  22)  into  the  present  solution  algorithm. 
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Figure  1  Flow  Field  Nomenclature 
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Figure  2  Flat-Plate  Triple-Deck  Pressure  Distributions 
(a)  Overall  Solution 
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Figure  2  Flat-Plate  Triple-Deck  Pressure  Distributions 
(b)  Trailing  Edge  Region 


81- 


-88-6 


- I  CO  (TRIPLE  DECK, 


Ax  =  0  2,  Ay -  04 


DISTANCE  FROM  TRAILING  EDGE  —  X 

Figure  3  Flat-Plate  Interacting  Boundary  Layer  Solutions 

(a)  Reynolds  Number  Effect  —  Forward/Reverse  Integration 
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Figure  3  Flat-Plate  Interacting  Boundary  Layer  Solutions 

(b)  Reynolds  Number  Effect  —  Forward  Integration 
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Figure  3  Flat-Plate  Interacting  Boundary  Layer  Solutions 
(c)  Global  Interaction  Effects 
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Figure  4  Finite  Reynoids  Number  Solutions  for  a  Flat-Plate  Trailing  Edge 
(a)  Surface  Pressure 
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Figure  4  Finite  Reynoids  Number  Solutions  for  a  Flat-Plate  Trailing  Edge 
(b)  Skin  Friction  and  Wake  Centerline  Velocity 
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Figure  4  Finite  Reynolds  Number  Solutions  tor  a  Flat-Plate  Trailing  Edge 
(c)  Displacement  Thickness 
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Figure  5  Laminar  Blunt  Trailing  Edge  Solutions 
(a)  Trailing  Edge  Geometry 
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Figure  5  Laminar  Blunt  Trailing  Edge  Solutions 
(b)  Displacement  Surface 
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Figure  5  Laminar  Blunt  Trailing  Edge  Solutions 

(c)  Skin  Friction  Parameter  and  Wake  Centerline  Velocity 
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Figure  5  Laminar  Blunt  Trailing  Edge  Solutions 
(e)  Pressure  Distribution 
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Figure  6  Reynolds  Number  Effects  for  Laminar  Blunt  Trailing  Edges 
(a)  Surface  Pressure 
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Figure  7  Turbulent  Flat  Plate  Results 
(a)  Pressure  Coefficient 
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Figure  7  Turbulent  Flat  Plate  Results 
(b)  Wake  Centerline  Velocity 
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Figure  8  Turbulent  Blunt  Trailing  Edge  Solutions 
(a)  Displacement  Surface 
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Figure  8  Turbulent  Blunt  Trailing  Edge  Solutions 
(b)  Pressure  Coefficient 
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Figure  9  Effect  of  Airfoil  Thickness  on  Turbulent  Trailing  Edge  Flow 
(a)  Displacement  Surface 
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Figure  9  Effect  of  Airfoil  Thickness  on  Turbulent  Trailing  Edge  Flow 
(b)  Skin  Friction  Parameter  and  Wake  Centerline  Velocity 
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Figure  9  Effect  of  Airfoil  Thickness  on  Turbulent  Trailing  Edge  Flow 
(c)  Pressure  Coefficient 
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